;***********************************************
; binning consistent with the "pdfxy" function
;
undef("pdfx_log")
function pdfx_log(x:numeric, nbin[1]:integer, opt:logical)
local nGood, nbins, xMin, xMax, mnmxint, xSpace \
     ,bin, pdf, nTot, XMIN, XMAX, nLoOut, nHiOut
begin

  if (opt .and. isatt(opt, "v510") .and. opt@v510) then
      pdf = pdfx_v510(x, nbin, opt) ; different binning
      return(pdf)
  end if

  if (nbin.le.2) then
      nbins = 25             ; default 
  else
      nbins = nbin
  end if

  nGood = num(.not.ismissing(x))
  if (nGood.lt.3) then
      print("pdfx: nGood="+nGood+" : Need more non-missing points")
      pdf = new( nbins, "double", getFillValue(x))
      return( pdf )
  end if

  XMIN = min(x)*1d0   ; min/max for ENTIRE array             
  XMAX = max(x)*1d0   ; force "double" [*1d0]

  xMin = 0.0d         ; prototype as double
  xMax = 0.0d

  if (opt .and. isatt(opt,"bin_min")) then  
      xMin = opt@bin_min                    ; user set
  else
      xMin = XMIN                           ; calculate
  end if

  if (opt .and. isatt(opt,"bin_max")) then  
      xMax = opt@bin_max                    ; user set
  else
      xMax = XMAX                           ; calculate
  end if

  if (opt .and. isatt(opt,"bin_nice")) then ; nice xMin, xMax
      outside = False
      if (isatt(opt,"bin_nice_outside")) then
          outside = opt@bin_nice_outside
      end if
      mnmxint = nice_mnmxintvl( XMIN, XMAX, nbins, outside)
      xMin    = mnmxint(0)
      xMax    = mnmxint(1)
      xSpace  = mnmxint(2)
    ;;nbins   = round( (xMax-xMin)/xSpace , 3) 
  end if

;;dbin        = (xMax-xMin)/(nbins-1)          ; 5.2.0   
  dbin        = (xMax-xMin)/nbins   
  binBound    = exp(fspan(log(xMin),log(xMax),nbins+1)) ;;xMin + ispan(0,nbins,1)*dbin
  binBound(nbins) = xMax                       ; avoid roundoff              
  binCenter   = (binBound(0:nbins-1) + binBound(1:nbins))*0.5d 

  binBoundMin = binBound(0)           
  binBoundMax = binBound(nbins)      

  pdf         = new( nbins, "double", getFillValue(x))
  pdf         = 0.0d

  do nb=0,nbins-1      
     pdf(nb) = num( x.ge.binBound(nb) .and. x.lt.binBound(nb+1) )
     if (nb.eq.(nbins-1)) then                 ; last bin 
         pdf(nb) = pdf(nb) + num( x.eq.binBound(nb+1) )   ; include last bound
     end if
  end do

  pdf!0    = "x"
  pdf&x    = binCenter
                                       ; max possible in data     
  nMax     = num(x.ge.XMIN .and. x.le.XMAX)
                                       ; actual number used
  nUse     = num(x.ge.binBoundMin .and. x.le.binBoundMax)

  nLoOut   = num(x.lt.binBoundMin)     ; number outliers
  nHiOut   = num(x.gt.binBoundMax)

  pdf      = 100d0*pdf/nMax            ; percent frequency

  pdf@bin_center     = binCenter
  pdf@bin_bounds     = binBound
  pdf@bin_bound_min  = binBoundMin
  pdf@bin_bound_max  = binBoundMax
  pdf@bin_spacing    = dbin            ; binBound(2)-binBound(1)
  pdf@nbins          = nbins
  pdf@nMax           = nMax
  pdf@nUse           = nUse
  if (nLoOut.gt.0 .or. nHiOut.gt.0) then
      pdf@nLoOut     = nLoOut       
      pdf@nHiOut     = nHiOut
  end if

  pdf@long_name      = "PDF"
  if (isatt(x,"long_name")) then
      pdf@long_name  = "PDF: "+x@long_name
  end if
  pdf@units          = "%"
  
  return( pdf )
end


